This R code downloads language data by census tract for Travis County, Texas, using 2013-2017 5-year ACS estimates. Then it maps those languages alongside Austin Public Health facilities. The goal is to reveal any areas that might be underserved.

Helpful websites:

library(tidyverse)
## ── Attaching packages ──────────────────────────────────────────────────────── tidyverse 1.2.1 ──
## ✔ ggplot2 3.1.0       ✔ purrr   0.3.2  
## ✔ tibble  2.1.1       ✔ dplyr   0.8.0.1
## ✔ tidyr   0.8.3       ✔ stringr 1.4.0  
## ✔ readr   1.3.1       ✔ forcats 0.4.0
## ── Conflicts ─────────────────────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
library(tidycensus)
library(viridis)
## Loading required package: viridisLite
library(sf)
## Linking to GEOS 3.6.1, GDAL 2.1.3, PROJ 4.9.3
library(readxl)
library(mapview)
options(tigris_use_cache = TRUE)

First read in Austin Public Health Locations data. This is found on the city open data portal (link?) but first you must break out lat and lon fields. Currently the dataset has lat and lon as part of the address field. Also added one more location, and then saved as an Excel file for the import.

Austin_Public_Health_Locations <- read_excel("Austin_Public_Health_Locations.xlsx")

Use tidycensus package to connect via Census API to download numbers of language speakers in each census tract. Include a variable for summary_var (total number of people in the census tract) and calculate percentages. Those percentage values will be used as the highlighted attribute for each census tract.

spanish <- get_acs(geography = "tract", variables = "C16001_003", year = 2017, state = 48, county = c("Travis"), summary_var = "C16001_001", geometry = TRUE, survey = "acs5") %>% mutate(pct = round(100 * (estimate / summary_est), digits = 1))
## Getting data from the 2013-2017 5-year ACS
chinese <- get_acs(geography = "tract", variables = "C16001_021", year = 2017, state = 48, county = c("Travis"), summary_var = "C16001_001", geometry = TRUE, survey = "acs5") %>% mutate(pct = round(100 * (estimate / summary_est), digits = 1))
## Getting data from the 2013-2017 5-year ACS
vietnamese <- get_acs(geography = "tract", variables = "C16001_024", year = 2017, state = 48, county = c("Travis"), summary_var = "C16001_001", geometry = TRUE, survey = "acs5") %>% mutate(pct = round(100 * (estimate / summary_est), digits = 1))
## Getting data from the 2013-2017 5-year ACS
otherasian <- get_acs(geography = "tract", variables = "C16001_030", year = 2017, state = 48, county = c("Travis"), summary_var = "C16001_001", geometry = TRUE, survey = "acs5") %>% mutate(pct = round(100 * (estimate / summary_est), digits = 1))
## Getting data from the 2013-2017 5-year ACS
frenchhaitian <- get_acs(geography = "tract", variables = "C16001_006", year = 2017, state = 48, county = c("Travis"), summary_var = "C16001_001", geometry = TRUE, survey = "acs5") %>% mutate(pct = round(100 * (estimate / summary_est), digits = 1))
## Getting data from the 2013-2017 5-year ACS
german <- get_acs(geography = "tract", variables = "C16001_009", year = 2017, state = 48, county = c("Travis"), summary_var = "C16001_001", geometry = TRUE, survey = "acs5") %>% mutate(pct = round(100 * (estimate / summary_est), digits = 1))
## Getting data from the 2013-2017 5-year ACS
russian <- get_acs(geography = "tract", variables = "C16001_012", year = 2017, state = 48, county = c("Travis"), summary_var = "C16001_001", geometry = TRUE, survey = "acs5") %>% mutate(pct = round(100 * (estimate / summary_est), digits = 1))
## Getting data from the 2013-2017 5-year ACS
otherindoeuropean <- get_acs(geography = "tract", variables = "C16001_015", year = 2017, state = 48, county = c("Travis"), summary_var = "C16001_001", geometry = TRUE, survey = "acs5") %>% mutate(pct = round(100 * (estimate / summary_est), digits = 1))
## Getting data from the 2013-2017 5-year ACS
korean <- get_acs(geography = "tract", variables = "C16001_018", year = 2017, state = 48, county = c("Travis"), summary_var = "C16001_001", geometry = TRUE, survey = "acs5") %>% mutate(pct = round(100 * (estimate / summary_est), digits = 1))
## Getting data from the 2013-2017 5-year ACS
tagalog <- get_acs(geography = "tract", variables = "C16001_027", year = 2017, state = 48, county = c("Travis"), summary_var = "C16001_001", geometry = TRUE, survey = "acs5") %>% mutate(pct = round(100 * (estimate / summary_est), digits = 1))
## Getting data from the 2013-2017 5-year ACS
arabic <- get_acs(geography = "tract", variables = "C16001_033", year = 2017, state = 48, county = c("Travis"), summary_var = "C16001_001", geometry = TRUE, survey = "acs5") %>% mutate(pct = round(100 * (estimate / summary_est), digits = 1))
## Getting data from the 2013-2017 5-year ACS
other <- get_acs(geography = "tract", variables = "C16001_036", year = 2017, state = 48, county = c("Travis"), summary_var = "C16001_001", geometry = TRUE, survey = "acs5") %>% mutate(pct = round(100 * (estimate / summary_est), digits = 1))
## Getting data from the 2013-2017 5-year ACS

This builds a new field in the data file, converting lat and lon fields into simple feature polygon geometry.

aph <- st_as_sf(Austin_Public_Health_Locations, coords = c("lon", "lat"),  crs = st_crs(spanish))

Produce a table of all languages and total # speakers.

totals <- tribble(~Lang, ~Total, "Spanish", sum(spanish$estimate), "Chinese", sum(chinese$estimate))
totals <- add_row(totals, Lang = "Vietnamese", Total = sum(vietnamese$estimate))
totals <- add_row(totals, Lang = "German", Total = sum(german$estimate))
totals <- add_row(totals, Lang = "French Haitian or Cajun", Total = sum(frenchhaitian$estimate))
totals <- add_row(totals, Lang = "Russian Polish Other Slavic", Total = sum(russian$estimate))
totals <- add_row(totals, Lang = "Other Indo-European", Total = sum(otherindoeuropean$estimate))
totals <- add_row(totals, Lang = "Korean", Total = sum(korean$estimate))
totals <- add_row(totals, Lang = "Tagalog", Total = sum(tagalog$estimate))
totals <- add_row(totals, Lang = "Arabic", Total = sum(arabic$estimate))
totals <- add_row(totals, Lang = "Other", Total = sum(other$estimate))
totals <- add_row(totals, Lang = "Other Asian", Total = sum(otherasian$estimate))

# Sort by # of speakers and print the totals.
totals <- totals %>% arrange(desc(Total))
totals
## # A tibble: 12 x 2
##    Lang                         Total
##    <chr>                        <dbl>
##  1 Spanish                     262579
##  2 Other Indo-European          20632
##  3 Chinese                      12395
##  4 Vietnamese                   11406
##  5 Other Asian                  11267
##  6 French Haitian or Cajun       5940
##  7 Korean                        4960
##  8 Arabic                        4866
##  9 Other                         4736
## 10 Russian Polish Other Slavic   3362
## 11 German                        3201
## 12 Tagalog                       1888

Remove tracts from the map with zero speakers. Need to expand this part. For now, just use a filter.

spanish <- filter(spanish, pct > 0)
otherindoeuropean <- filter(otherindoeuropean, pct > 0)
chinese <- filter(chinese, pct > 0)
vietnamese <- filter(vietnamese, pct > 0)
otherasian <- filter(otherasian, pct > 0)
arabic <- filter(arabic, pct > 0)
frenchhaitian <- filter(frenchhaitian, pct > 0)
korean <- filter(korean, pct > 0)
other <- filter(other, pct > 0)
russian <- filter(russian, pct > 0)
german <- filter(german, pct > 0)
tagalog <- filter(tagalog, pct > 0)

Set the color brewer ramp colors. We pick the six monochromatic ramps in RColorBrewer.

library(RColorBrewer)
reds.pal <- colorRampPalette(brewer.pal(9, "Reds"))
purples.pal <- colorRampPalette(brewer.pal(9, "Purples"))
oranges.pal <- colorRampPalette(brewer.pal(9, "Oranges"))
greys.pal <- colorRampPalette(brewer.pal(9, "Greys"))
greens.pal <- colorRampPalette(brewer.pal(9, "Greens"))
blues.pal <- colorRampPalette(brewer.pal(9, "Blues"))

Build the first map, with the top five languages in the city. Map is shown below; larger version on a standalone html page can be found here.

map <- aph %>% filter(`Facility Type` == "WIC Clinic") %>% mapview(zcol = "Facility Name", legend = FALSE, layer.name = "WIC Clinics", col.regions = "green", homebutton = FALSE)
map <- map + aph %>% filter(`Facility Type` == "Neighborhood Center") %>% mapview(zcol = "Facility Name", legend = FALSE, layer.name = "Neighborhood Centers", col.regions = "blue", homebutton = FALSE)
map <- map + aph %>% filter(`Facility Type` == "Other Location") %>% mapview(zcol = "Facility Name", legend = FALSE, layer.name = "Other Locations", col.regions = "yellow", homebutton = FALSE)
map1 <- map + mapview(spanish, zcol = "pct", col.regions=reds.pal, homebutton = FALSE)
map1 <- map1 + mapview(otherindoeuropean, zcol = "pct", col.regions=purples.pal, homebutton = FALSE)
map1 <- map1 + mapview(chinese, zcol = "pct", col.regions=greens.pal, homebutton = FALSE)
map1 <- map1 + mapview(vietnamese, zcol = "pct", col.regions=blues.pal, homebutton = FALSE)
map1 <- map1 + mapview(otherasian, zcol = "pct", col.regions=oranges.pal, homebutton = FALSE)
map1 <- map1 + mapview(arabic, zcol = "pct", col.regions=purples.pal, homebutton = FALSE)

mapshot(map1, url = paste0(getwd(), "/docs/map1.html"))
map1

Build the second map, of the bottom seven languages. Map below; final map can also be found here.

map2 <- map + mapview(frenchhaitian, zcol = "pct", col.regions=reds.pal, homebutton = FALSE)
map2 <- map2 + mapview(korean, zcol = "pct", col.regions=purples.pal, homebutton = FALSE)
map2 <- map2 + mapview(arabic, zcol = "pct", col.regions=oranges.pal, homebutton = FALSE)
map2 <- map2 + mapview(other, zcol = "pct", col.regions=greens.pal, homebutton = FALSE)
map2 <- map2 + mapview(russian, zcol = "pct", col.regions=blues.pal, homebutton = FALSE)
map2 <- map2 + mapview(german, zcol = "pct", col.regions=reds.pal, homebutton = FALSE)
map2 <- map2 + mapview(tagalog, zcol = "pct", col.regions=purples.pal, homebutton = FALSE)

mapshot(map2, url = paste0(getwd(), "/docs/map2.html"))
map2